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2 ■ Abstract 

, Context. Herbig-Haro objects are characterized by a complex knotty morphology detected mainly along the axis of protostellar jets in a wide 
1-1 range of bands: from radio to IR to optical bands, with X-rays knots also detected in the last few years. Evidence of interactions between knots 
, ' formed in different epochs have been found, suggesting that jets may result from the ejection of plasma blobs from the stellar source. 

Aims. We aim at investigating the physical mechanism leading to the irregular knotty structure observed in protostellar jets in different 
wavelength bands and the complex interactions occurring among blobs of plasma ejected from the stellar source. 

Methods. We perform 2D axisymmetric hydrodynamic numerical simulations of a randomly ejected pulsed jet. The jet consists of a train of 
, blobs which ram with supersonic speed into the ambient medium. The initial random velocity of each blob follows an exponential distribution. 
£S| ■ We explore the ejection rate parameter in order to derive constraints on the physical properties of protostellar jets by comparison of model 
results with observations. Our model takes into account the effects of radiative losses and thermal conduction. 

Results. We find that the mutual interactions of blobs ejected at different epochs and with different speed lead to a variety of plasma 
components not described by current models of jets. The main features characterizing the random pulsed jet scenario are: single high speed 
knots, showing a measurable proper motion in nice agreement with optical and X-rays observations; irregular chains of knots aligned along 
r* ■ the jet axis and possibly interacting with each other; reverse shocks interacting with outgoing knots; oblique shock patterns produced by 
# the reflection of shocks at the cocoon surrounding the jet. All these structures concur to determine the morphology of the jet in different 
k> wavelength bands. We also find that the thermal conduction plays a crucial role in damping out hydrodynamic instabilities that would develop 
^ ■ within the cocoon and that contribute to the jet breaking. 
| 
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1. Introduction 

Herbig-Haro (HH) objects characterize the emission regions 
near protostars where the accretion process is still ongoing, of- 
ten allowing us to identify the presence of a young stellar object 
(YSO) not observable otherwise. The HH objects are observed 
in all wavelength bands (from radio and IR to optical and UV, 
and a lso in the X-rays; see the review of iReipurth & Ballyl 
2001) as a chain of knots within the supersonic protostellar jets. 
The knotty structure has been frequently observed in the optical 
band in many HH objects and has been discovered for the first 
time also in X-rays in HH 154 dFavata et al .1120061 iBonito et al.1 
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2008b . one among the nearest and most luminous X-ray emit- 
ting protostellar jets. The HH 154 jet is one of the best studied, 
being observed through multi- wavelength and multi-epoch ob- 
servations with several instruments, among them the Hubble 
Space Telescope (HST; iFridlund & Li seau 1998; Bonito etaD 
120081) . Chandra/AC IS-I dBallv et al.ll2003l: Ir^vata et alJl2006h. 
and XMM-Newton dFavata et al.ll2002h . 

The origin of chains of optical knots commonly observed 
in jets have been traditionally interpreted through models of 
jet eje cted with a sinusoidal variable velocity (e.g. Rag a et al 



2007k . In fact, in these models, the internal shocks resulting 



from plasma flowing with different velocity lead to a nodular 
structure, constituted by a regular chain of knots with individ- 
ual proper motions. An alternative scenario has been proposed 
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bv lYiraketaDd2009h who model, through 3D simulations, the 
appearing of clumps within the jet beam (with random loca- 
tions, sizes, and speeds), and explore the possibility of intrin- 
sic density heterogeneity of the jet. Original models of jets ac- 
counting for stochastic jets and random tim e variability have 
been previously discus sed by Roberts ( 1986|) (in the context of 
extragalactic jets) and iRagal (I1992|) (in the limit of larg e dis- 
tances from the stellar source). Also lVitorino et al.1 (I2002|) sug- 
gested a disk/jet launching model with a combination of a con- 
tinuous and a random component of the velocity, thus justifying 
the random speed models. 

Recently, multi- wavelength observations of HH 154 col- 
lected in 2005 with both Chandra and HST showed evidence 
of mutual interactions between bright knots, for med proba 



bly at different epochs and with different speeds (IBonito et al 



2008) . These features s eem to be common also in other ob 
jects dGaratti et al . 2009) and are not predicted in mo dels of jets 
with s inusoidal velocity. Moreover as discussed by lYirak et al 



(2009), the observed knots could be interpreted erroneously as 
the product of a sinusoidal variable velocity jet. 

The aim of this work is to investigate, through detailed 
numerical modeling, the origin and dynamics of the irregular 
knotty structure observed in several bands (and, in particular, 
in X-rays) in protostellar jets and the possible interactions be- 
tween knots. Moreover, we aim at deriving general model pre- 
dictions to be compared with observations of X-ray emitting 
protostellar jets (DG Tau, HH 2, etc.) to understand the physi- 
cal mechanisms leading to X-ray emission in different objects. 

In our previous studies, we have investigated the origin 
of X-ray emission from protostellar jets , through a model 
of a continuously ejected supersonic jet dBonito et al.1 12004 



IBonito et al.ll2007[) . Such a model predicts the propagation of a 
compact X-ray source through the circumstellar medium with 
spectral behavior, X-ray luminosi ty, and proper mot ion in good 
agreement with the observations dFavata et al.ll2002|) . However, 
the model does not reproduce the complex knotty morphol- 
ogy and variability of the X-ray source detected in HH 154 
dFavata et al.ll2006|) . 

Here we suggest an improvement to our previous model 
by considering the scenario of a pulsed jet interacting with an 
initially homogeneous ambient medium. The model presented 
here describes a jet constituted by a train of plasma blobs ran- 
domly ejected by the stellar source. The initial velocity of our 
model follows an exponential distribution and this choice cor- 
responds to a distribution of many low-speed blobs and few 
blobs with high initial ejection velocity. The blobs are expected 
to interact with each other leading to the formation of knots 
with individual proper motions, possibly observable in various 
bands. In the following, we will refer to blobs to consider the 
plasma ejected from the stellar source and to knots to consider 
the observable structures within the jet. 

We plan to present our results in two papers. In the present 
paper we describe the dynamics and energetic of the pulsed 
jet/ambient interaction, exploring the influence of the ejection 
rate of plasma blobs. In the second paper (Bonito et al. 2009, in 
preparation) we will focus on the variability and morphology of 
X-ray emitting features predicted by the model as they would 
be observed with Chandra/ ACIS -I and with XMM/Newton. 



The paper is organized as follows: in Sect. [2] we describe 
the numerical model of the pulsed jet, its initial setup, and the 
exploration of the parameter space; in Sect. [3] we analyze the 
hydrodynamic evolution of the model, focusing on the internal 
structures occurring within the jet, and the physical effects due 
to the thermal conduction process; in Sect. [4] we discuss the 
results and their physical interpretation, and in Sect. [5] we draw 
our conclusions. 

2. The model 

The model describes a protostellar jet ramming with a super- 
sonic variable ejection speed into an initially homogeneous am- 
bient medium. The fluid is assumed to be fully ionizecQ with a 
ratio of specific heats 7 = 5/3. The model takes into account 
the radiative cooling from optically thin plasma, and the ther- 
mal conduction (including the effects of heat flux saturation). 
The propagation of the jet through the surrounding medium is 
modeled by solving numerically the time-dependent hydrody- 
namic equations: 



dp 

dp\ 
~di 

dpE 
~df 

with 



V • p\ = 
f- V • p\\ + Vp = 

+ V • (pE + p)y = - V • q - n e n H P(T) 



(1) 

(2) 
(3) 



p = (7 - l)pe , 



E = e 



where p is the pressure, E the total gas specific energy (internal 
energy, e, and kinetic energy) respectively, p is the mass den- 
sity, t the time, v the plasma velocity, q the heat flux, n e and 
Tin are the electron and hydrogen density, respectively, P(T) 
is the radiative losses function per unit emission measure^ T 
the plasma te mperature. 



Following lDalton & Balbusl (11993). we use an interpolation 
expression for the thermal conductive flux which al lows for a 
smoot h transition between the Spitzer ( Spitzej[l962 ) and satu- 
rated ( Cowie & McKeetfl977l) conduction regimes 



q 

where 



(-+- 

\ ^spi <2sat 



g 8 pi = -k(T)VT , 

where k(T) = 9.2 x 10~ 7 T 5 / 2 erg s" 
mal conductivity, and 

<2sat = -sign(VT)b(j)pcl , 



K- 



cm 



(4) 



(5) 

1 is the ther- 



(6) 



1 See lBonito et al.ll2007l . 

2 P(T) is described by a functional form that takes into account 
free-free, bound-free, bound-bound and 2 photons emission (s ee 
iRavmond & Smith 1977; Mewe et al. 1985; Kaastra & Mewe 2000). 
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where 0-0.3 dGiulianilll984]:lBorkowski et al.ll 19891 and ref- 
erences therein) and c s is the isothermal sound speed. 

T he hydrodynamic e quations are solved using the FLASH 
code ( Fryxell et al.l l2000) with customized numerical modules 



for the op tically thin radiativ e losses and the thermal conduc- 
tion (see lOrlando et al] [20051 for the details of the implemen- 



tation of the thermal conduction). The FLASH code uses a di- 
rectionally split Piecewise Parabolic Meth od (PPM) solver to 
handl e compressible flows with shocks (IColella & Woodward 
1984 . the PARAM ESH library to handle adaptive mesh re- 
finement (iMacNeice et al.1 12000b . and the Message-Passing 
Interface library to achieve parallelization. 

2.1. Numerical setup 

We solve the hydrodynamic equations in 2D, using cylindrical 
coordinates in the plane (r, z), and assuming axisymmetry with 
the jet axis coincident with the axis of symmetry (z). The com- 
putational grid size is (2000 AU x 6000 AU). The maximum 
spatial resolution achieved by our simulations is « 8 AU, as 
determined with the algorithm of the PARAMESH library for 
4 refinement levels, which corresponds to covering the initial 
jet radius with 4 points at the maximum resolution [j. 

We follow the evolution of the jet for a time interval rang- 
ing between 100 and 400 yrs, depending on the set of model pa- 
rameters. The pulsed jet consists of a train of blobs, each lasting 
for 0.5 yr, with an ejection rate corresponding to a time interval 
between two consecutive blobs At. Each blob is ejected with a 
random velocity directed along the jet axis following an expo- 
nential distribution 



v(r = 0, t) = VjMAx(b — log(rand(t) — a)) 
where 

1 



b = log 



1 - e~ 



(7) 



(8) 
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Figure 1. Initial jet velocity of each blob as a function of time 
for the At = 2 yr case (run LJ2). 



where v(r = 0, t) is the on-axis velocity, v is the ambient- to-jet 
density ratio, rj is the initial jet radius and w = 4 is the steep- 
ness parameter for the shear layer which allows a smooth tran- 
sition of the kinetic energy at the interface between the jet and 
the ambient medium. This lay er avoids the grow th of random 
perturbation at the jet border dBodo et al.lll994|) . The density 
profile as a function of the radial direction is 



p(r) = pAv 



cosh(r/rj) 1 



(11) 



where p- } is the jet density (see also the velocity and density 
profiles as a func tion of the radial distance from the axis in 
Bonito et al.lE007[ Fig. 1). 

Reflection boundary conditions are imposed along the jet 
axis (consistent with the adopted symmetry), inflow boundary 
conditions at the base for r < r 3 ; , where r is the radial distance 
in cylindrical coordinates and Tj is the jet radius, and outflow 
boundary conditions elsewhere. 



1 ' 



(9) 2.2. Parameter space 



and rand is a number ranging between and 1 obtained using a 
random number generator. The maximum velocity of the blobs 
is VjMAX = 4680 km/sec and corresponds to the highly super- 
sonic initial Mach number M = 1000. Fig. Q] shows the initial 
velocity values of each blob ejected with an ejection rate cor- 
responding to At = 2 yr. We generated the random seed once 
for all the runs so the values follow always the same sequence, 
only the rate is different. 

The initial random velocity of each blob is along the z axis, 
coincident with the jet axis, and its radial profile is 



v(r,t) 



v(r 0,t) 



v cosh(r /rj) w — {y — 1) 



(10) 



3 The maximum resolution corresponds to ~ 8 times the resolution 
of Chandra/ ACIS -I observations and to about twice the spatial resolu- 
tion of HST, 0.4"^ 60 AU at 150 pc and 0.1"« 15 AU at 150 pc, i.e. 
at the distance of the nearest star forming region in Taurus. 



The model parameters are derived from t he optical and X- 
ray data analysis of protostellar jets (e.g. iPodio et al. 2006; 



Bonito et al. 2008: Favata et a I. 2006) and from the results of 
previous models "dBonito et alJl2004l2007h . Tab Q] summarizes 
the physical parameters characterizing the simulations. 

The exploration of the parameter space mainly focuses on 
the ejection rate of the blobs composing the jet. In all our sim- 
ulations, we impose the initial radius and density of the blobs, 
namely rj « 30 AU and rij = 500 cm -3 , respectively. For 
the sake of completeness, we explore both the initial light jet 
scenario (a jet initially less dense than the unperturbed ambi- 
ent medium), and the initial heavy jet scenario (a jet initially 
denser than the ambient medium); we consider, therefore, ini- 
tial ambient- to-jet density ratio v = n a /rij = 0.1, 10 in the 
heavy and light jet cases, respectively. The initial blob tem- 
perature Tj ranges between 10 2 and 10 4 K, while the ambient 
temperature, T a , is derived assuming initial pressure equilib- 
rium between the jet and the unperturbed ambient. 
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Table 1. Summary of the initial physical parameters character- 
izing the simulations: At is the time between two consecutive 
blobs, v is the ambient- to-jet density contrast, M is the ini- 
tial jet Mach number, v- } is the velocity of each ejected blob 
(here we indicate the range of values randomly generated in 
our model), Tj is the initial jet temperature, n a is the ambient 
density. In all the models we assume: rj = 30 AU as initial jet 
radius, rij = 500 cm -3 as initial jet density, and T a = 10 3 K 
as initial ambient temperature. 



Model 


At 


V 


M 






n a 




yr 






[kms -1 ] 


[10 4 K] 


[cm- 3 ] 


LJ0.5 


0.5 


10 


1000 


10 - 4700 


10 4 


5000 


LJ2 


2 


10 


1000 


10 - 4700 


10 4 


5000 


LJ8 


8 


10 


1000 


10 - 4700 


10 4 


5000 


LJ2RD 


2 


10 


1000 


10 - 4700 


10 4 


5000 


HJ2 


2 


0.1 


300 


10 - 1400 


10 2 


50 



The initial jet Mach number values are M = Vj/c a = 300 
in the heavy jet case and 1000 in the light jet case, and corre- 
sponds to a maximum velocity of the first blob of about 5000 
km/s. Note that we are interested on the evolution of the pulsed 
jet after a transient phase, i.e. after the first ejected blob of 
plasma has perturbed the ambient (which is no more homo- 
geneous). We require, therefore, a high speed of the first blob 
in order to perturb quickly the computational domain (6000 
AU). The speed of each ejected blob ranges between about 10 
km/s and 5000 km/s during the evolution of the pulsed jet. It 
is worth noting that the apparently high values of the maxi- 
mum ejection velocity cannot be compared directly to the ob- 
servable speeds of the knots within the protostellar jets. In fact 
we have chosen the initial speeds which lead, once the interac- 
tion between the ejected plasma and the ambient medium has 
taken place, to the observed knots speed of a few hundred km/s 
dFridlund & Liseaull 19981: [Bonito et alJl2008h . 

Since the optical knots in HH jets, in general, show time- 
scale variability of a few years (see morph ological evolution o f 
the HH 154 jet in about 2 yr discussed in lBonito et"aDl2Q08h . 
with new knots appearing at the base of the jet in a bout 2 yr 
(e.g. HH 1 and HH 34, see the movies discussed in iHartiganl 
2003), we have chosen an initial ejection rate parameter corre- 
sponding to a time interval At = 2 yr between two consecutive 
blobs. Then we explore this parameter taking into account a 1/4 
value (At = 0.5 yr) and a x4 value (At = 8 yr) in order to con- 
strain the ejection rate from the comparison between our model 
predictions and the observations. 

3. Results 

3.1. Hydrodynamic evolution 

In all our simulations, the jet propagates through the homoge- 
neous ambient medium, forming a cocoon, and shows a com- 
plex density structure inside. The ejection rate determines the 
energy powering the jet and, in fact, the maximum distance 
from the driving source reached at a given epoch in runs LJ0.5, 
LJ2, and LJ8 strongly varies. As shown in Fig. [2 the low ejec- 
tion rate case (LJ8) travels through 3000 AU in 100 yr, while 
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Figure 2. 2D spatial distributions of the density of the 
jet/ambient interaction in runs LJ2 (left semi-panel) and LJ8 
(right semi-panel). The head of the jet leaves the domain after 
about 90 yr in LJ2 whereas in LJ8 it does not reach the end of 
the domain (6000 AU) in 100 yr. 



the high ejection rate cases reach the 6000 AU distance in less 
than 90 yr (LJ2) or even less than 50 yr (LJ0.5), being the ki- 
netic power released in these cases much higher than in LJ8. 

Although the ambient medium is initially homogeneous, 
it becomes quickly inhomogeneous as the first and consecu- 
tive blobs ejected by the stellar source interact with it. Fig. [3] 
shows a cut along the jet axis of the density (upper panel) and 
of the temperature (lower panel) distributions. The whole com- 
putational domain (6000 AU) is filled with several structures, 
with the density varying by more than two orders of magni- 
tude. Also the thermal conditions vary along the jet axis, so 
that the sound speed and the Mach number are not uniform. As 
a result, strong variations of the shock velocities (v 3 h oc p -1 / 2 ) 
and of the post-shock properties are determined due to the vary- 
ing pre- shock conditions. Note that, at variance with the model 
of a continuous jet, in our pulsed jet model the density con- 
trast between the ambient medium and the blob can vary during 
the jet/ambient evolution. As discussed above, in a few years, 
the medium is so complex that, even if the initial density of 
the ejected blobs does not vary (nj = 500 cm -3 , see the red 
dashed line superimposed in the upper panel of Fig. 0, the am- 
bient density strongly varies during the evolution thus leading 
to several different scenarios in which the ejected blob could be 
denser or less dense than the surrounding medium (see upper 
panel in Fig.0. 

It is worth emphasizing that the random velocity of the 
blobs ejected at different epochs together with the strongly in- 
homogeneous medium in which each blob (except the first one) 
propagates, lead to complex mutual interactions of the blobs. 
The frequency of mutual interactions within the jet depends on 
the ejection rate, being the largest in run LJ0.5 and the lowest 
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Figure 3. Cut along the jet axis of the spatial distributions of 
density (upper panel) and of temperature (lower panel) in run 
LJ0.5: the chain of knots within the jet due to previously ejected 
blobs lead to a complex structure of the medium where the 
new blobs are ejected. The red dashed line superimposed in the 
density profile indicates the initial jet density value (n- } — 500 
cm -3 ). 



in run LJ8. These interactions result in a rich variety of plasma 
structures that cannot be described by current models of jet and 
that determine an irregular pattern of knots. In the following, 
we illustrate the common features that may contribute in de- 
termining the morphology of the pulsed jet by analyzing the 
density and temperature spatial distributions and evolution. 

The most evident feature predicted by our model is an irreg- 
ular chain of knots aligned along the jet axis. Depending on the 
model parameters, the knots may emit mainly either in the op- 
tical band or in the X-ray band. In the former case, the chain of 
(optical) knots can be highlighted by deriving density maps of 
plasma with temperature ranging between (5000 — 100000) K, 
which is a proxy of optical emission. To this end, we first recon- 
struct the 3D spatial distributions of mass density and tempera- 
ture, according to the symmetry of the problem; then, assuming 
the jet traveling perpendicularly to the line-of-sight (LoS), we 
integrate the density along the LoS and derive the 2D density 
maps. As an example, Fig. [4] shows an enlargement of the base 
of the jet for run HJ2, i.e. in the heavy jet scenario (see Tab. 



[B, where a chain of (optical emitting) knots within the jet is 
well visible. These knots show a detectable proper motion in a 
few years, in good agreement with typical observations of the 
moving optical knots within protostellar jets. In particular, by 
measuring the spatial displacement of the center of the circu- 
lar region identifying the E knot in Fig. 01 in the time interval 
(45 — 50) yr, we derive an average knot speed along the jet axis 
of about 600 km/s, in nice agreement with the values obtained 
from th e optical observations of the knots within protostel- 
lar jets dFridlund et al.ll2005l: lionito et al.lE008l lEisloffel et~aT 



1 1994b . In this case, we find evidence of knotty structures anal- 
ogous to those observed in protostellar jets, with a chain of 
several subsequent knots in a few hundreds AU from the stellar 
source. 

Similar patterns of knots but emitting mainly in the X-ray 
band (see Bonito et al., 2009, in preparation) are found in the 
light jet scenario when the ejection rate is high (At = 0.5, 2 
yr), as shown in the 2D spatial distribution of mass density in 
Fig. [5j In particular several dense structures are well visible 
within the jet, mainly on the jet axis. 

The spatial separation between two consecutive knots is 
neither uniform in space nor constant in time in all the differ- 
ent cases analyzed because of the variation of the blob speed 
and of the medium surrounding the blob. This leads to a com- 
plex pattern consistent with obs ervations of several protostellar 
jets, among them HH 154 (see lBonito et al . 2008). Sinusoidal 
ejection velocity models usu ally show a c hain of knots with an 
almost regular pattern (IRaga et al .1120071) within the jet which 
does n ot reproduce the obse rved struct ure of complex jets , (e.g. 
HH l. lReipurthetaDl2000l: HH lll lHartigan eTalhoOll) . 

The knots resulting from the interactions among ejected 
blobs have, in general, different proper motions and may ulti- 
mately interact with each other. In some cases, high speed knots 
can interact with the head of the jet, possibly overtaking the co- 
coon, and thus traveling into the unperturbed ambient medium. 
A random variability of the injection velocity of the blobs and 
the consequent catching -up processes between knots have been 
discussed previously by lRagal (11992b . This effect is evident in 
Fig. |4l where the initially separated C and B knots merge in a 
single structure (named C + B in the figure). Analogously the 
density cut along the jet axis indicate that the two peaks C and 
B at t « 45 yr of evolutionary stage (upper panel of Fig. [6]), 
correspond to a single peak C + B after 5 yr (lower panel of 
Fig. [6]). These features have been obs erved in several HH jets 
dBonito et al.ll2008l:lGaratti et al.ll2009h . 

Analogous clump-clump interactions and high speed knots 
overtaking the head of the jet (includi ng also the c ase o f off- 
axis structures) have been found by lYirak et all d2009h . al- 
though these authors have not considered any proxy of opti- 
cal emissio n as d one here. Note, however, that the knots in 
Yirak et al. I d2009h are the result of imposed inhomogeneities 



of the jet flow (describing the effect of possible instabilities) 
and their characteristic size is lower than the jet diameter (sub- 
radius), whereas in our model the knots are the result of the 
pulsed ejection of the jet (possibly related to accretion phenom- 
ena) andjhdr^zej^comparable to the jet diameter. Note, also, 
that lYiraket all (12009) considered ranges of temperatures and 
velocities smaller than those discussed here and their model 
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Figure 4. Density maps of plasma with temperature ranging between (5000 — 100000) K (a proxy of optical emission) for the 
model HJ2. The green circles (of radius « 80 AU) superimposed to each frame mark the position of a specific knot as it evolves 
in time: the proper motion is detectable in a few years, the time-lapse between two consecutive frames is 1 yr. 
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Figure 5. 2D spatial distribution of the density of the 
jet/ambient interaction in model LJ0.5: a chain of several knots 
arises due to blobs ejected at different epochs. 



does not take into account the thermal conduction effects. As 
discussed below in Sect. l3.2l our simulations show that the ther- 
mal conduction plays a crucial role in damping out the hydro- 
dynamic instabilities that would develop within the jet, leading 
to a faster and more collimated jet. 

Another important feature produced by our model is associ- 
ated with the formation of reverse shocks moving in the oppo- 



site direction with respect to the ejected blobs and interacting 
with these (see Fig. [7]). This interaction produces, under cer- 
tain circumstances, the formation of shocks with no detectable 
proper motion over time scales of a few ye ars, as observed, fo r 
instance, in the X-ray source of HH 154 i Favata et~aT1l2006[) . 
The collision between knots and reverse shocks is a feature 
common to all the ejection rates discussed here and occurs sev- 
eral times during the time covered in our simulations. As an 
example, Fig. [7] shows one of these occurrence for run LJ0.5. 
In the first panel, the reverse shock is located « 1200 AU above 
the source driving the jet and the knot traveling towards it is ini- 
tially at about 200 AU from the source. As the jet evolves, the 
knot reaches a distance of ~ 800 AU from the source, with the 
reverse shock just « 200 AU above it. After about two years 
the two structures interact at a distance of « 1000 AU from 
the source and the shock location does not change appreciably 
over a timescale of few years. 

A more complex morphology arises when the reverse shock 
is reflected obliquely by the cocoon, leading, in this case, to an 
oblique shock pattern. Fig. [8] shows the evolution of an oblique 
shock in a few years, moving from the external cocoon toward 
the jet axis, as indicated by the arrows superimposed on each 
frame. Note that the propagation of the oblique shock causes 
the formation of a dense and well collimated core within the jet 
itself as evident in the last frame of Fig.[8j 

3.2. Effects of the thermal conduction 

Our model takes into account both the radiative losses from an 
optically thin plasma and the thermal conduction, in Spitzer or 
in saturated regimes. It is interesting, therefore, to investigate 
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Figure 6. Cuts along the jet axis of the first (upper panel) and 
the last (lower panel) frame of Fig. HI corresponding to 45 yr 
and 50 yr since the beginning of the jet/ambient interaction. In 
a few years, the proper motion of the knots within the jet, the 
interaction between knots, and the generation of new knots are 
well visible. 



the effects of the thermal conduction in the dynamics and en- 
ergetic of the jet. These effects are highlighted by comparing 
simulations either with or without the thermal conduction. To 
this end, we performed an additional run with setup identical 
to that used for run LJ2 but without the thermal conduction 
(run LJ2RD). We find that, in run LJ2RD, the maximum dis- 
tance from the stellar source is significantly smaller than that 
in run LJ2. As an example, Fig. [9] shows the temperature dis- 
tributions in the two cases after the first blob has reached the 
upper boundary of the computational domain, namely at t ~ 90 
yr in LJ2 and at t « 350 yr in LJ2RD. We found, therefore, that 
the jet travels more slowly into the ambient medium if the ther- 
mal conduction effects are inhibited. In fact, in LJ2RD, the jet 
evolution is dominated by thermal and hydrodynamic instabil- 
ities (see right panel in Fig. [9]) that are known to contribute to 
the jet breaking. Conversely, these instabilities are damped out 
by the thermal conduction in run LJ2, and the jet propagates to 
larger distances. 

The comparison between the two temperature maps of runs 
LJ2 and LJ2RD in Fig. [9] indicate also significant differences 



in the morphology due to the presence of the thermal conduc- 
tion. In fact, beside the different evolutionary stages needed in 
the two cases in order to reach the top of the computational 
domain, it is also evident how the thermal conduction leads 
to a smooth and more collimated spatial distribution of tem- 
perature. The radiative case LJ2RD shows a significant radial 
spread of the temperature distribution, accounting for the de- 
celerated material traveling into the instabilities. 

4. Discussion 

4. 1. Observable proper motion 



As discussed in iBonito et al.l (I2007|) . the high initial jet ve- 
locity is not an observable since we can measure just the 
speed of formed knots, i.e. the velocity of the plasma struc- 
tures formed by the interactions among the blobs or of the 
blobs with the inhomogeneous surrounding medium. In other 
words, we can observe the effect of jet/ambient interaction, 
not the launching region itself. We are interested, therefore, in 
deriving the average speed of the observable knots, determin- 
ing the internal density structure of the jet. To this end, Fig. 
TO] shows the density- weighted average of the z— component 
of the plasma bulk velocity within a distance of 10 pixel^fl 
from the jet axis for the three ejection rates analyzed (runs 
LJ0.5, L2, L8). The velocity values are consistent with typical 
knots speeds observed within optical jets ( few hundreds km/s , 
Eisloffel&Mu^dtlll998l: iBallv et alJEool IBonito et alJl2008h . 
In particular we find a very good agreement for the At = 2 yr 
case and the speed value decreases with increasing At param- 
eter. 

Following our model results, we can make predictions 
about the initial ejection velocity leading to the observed speed. 
In order to verify if our assumed values are reasonable, we can 
compare the mechanical luminosity derived from our simula- 
tions with the observed bolometric luminosity of the source of 
the HH 154 jet. In particular, from the values of the initial ve- 
locity used here and show n in Tab. [T] and following the relation 
(18) in lBonito et al.1 d2007 ). the kinetic power derived from our 
model ranges between L me ch ~ 10 -7 — 10 L , i.e. always 
below the observed bolometric lu minosity of the so urce of the 



HH 154 jet, L boi « 40 L© (see IBonito etaD 120071 Tab. 1). 



Thus we conclude that the initial blob velocity values assumed 
in our models lead to reasonably small kinetic power. 



Finally, following the relation (17) from IBonito et al 



(2007), we derive a jet mass loss rate ranging between M ~ 
10~ n — 10~ 9 M^ lyv. Once again, as discussed in details in 
Bonito etaD (2007), these results concerning both the kinetic 
power and the mass loss rate are several orders of magnitude 
lower than those observed in CO outflows and it seems reason- 
able that the simulated protostellar jets cannot drive molecular 
outflow. 



4.2. Ejection rate of the blobs 

Our simulations show that observable knots are mainly pro- 
duced by the mutual interactions among blobs of plasma. We 



Corresponding to a physical radius of « 70 AU. 
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Figure 7. Evolution in a few years of the 2D spatial distributions of the density of the jet/ambient interaction in model LJ0.5: 
reverse shocks (highlighted by the cyan circle in each panel) form and travel in the direction opposite to the ejected blobs 
(highlighted by the magenta circle in each panel) and interact with them (white circle superimposed on the last panel). 
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Figure 8. Model LJ8: evolution of the 2D spatial distributions of the density of the jet/ambient interaction (the first frame is after 
70 yr and the last after 75 yr since the beginning of the jet/ambient interaction). Oblique shocks form (marked by the arrows in 
each panel) as the ejected blob interacts with the lateral cocoon. 



expect, therefore, that the frequency of interactions between 
knots depends on the ejection rate of the blobs (see Sect. [2]). In 
fact, from our exploration of the parameter At, we find that the 
knot interactions occur quite frequently in runs with high ejec- 
tion rate (LJ05, and LJ2), whereas it is almost absent in the low 



ejection rate case (LJ8) over the timescale analyzed here. This 
result can be understood from the following simple estimation. 



Bonito et al.l (I2007|) showed that the observable velocity of 
a knot is a fraction of the initial jet velocity. This result is con- 
firmed by our simulations: as discussed in Sect. I4.ll Fig. [TO] 
shows that the maximum average speed of the knots along the 
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Figure 9. 2D spatial distributions of the temperature of the 
jet/ambient interaction in the At = 2 yr case: with radiative 
losses and thermal conduction after about 90 yr (left panel) and 
the pure radiative case after about 350 yr (right panel) of evo- 
lution. Ignoring the effects due to the thermal conduction, the 
whole computational domain is fully perturbed along the jet 
axis after more than 350 yr. 
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Figure 10. Average speed along the jet axis for the three runs 
LJ0.5, LJ2, and LJ8: the At = 0.5 yr case is in black, the 
At = 2 yr case is in red, and the At = 8 yr case is in green. The 
values of the knots speed derived from the runs are consistent 
with optical observations. 

jet axis is a fraction of the initial average blob velocity 2000 
km s _1 ). The time at which we expect to observe a collision 
between two consecutive knots is 

tco\\ = knot (At)knot (12) 
A^knot 

where A^knot and (At) knot are the average relative speed and 
the time interval between two consecutive knots, respectively. 
Both ^knot and A^knot are calculated assuming the same de- 
celeration factor (since each knot velocity is a fraction of the 
initial blob speed); the latter is evaluated by considering only 
the blobs ejected with higher speed with respect to the previous 



ones. Since knots are produced mainly by the mutual interac- 
tion between blobs and between blobs and the ambient, we ex- 
pect that (At) knot > At. Considering At instead of (At)k no t 
in Eq.[l2l we have, therefore, a lower limit on the collision time 
between knots. We derive that the interactions rises on average 
over about t co ii = 10, 30, and more than 100 yr in runs LJ0.5, 
LJ2, and LJ8, respectively. This estimate explains the presence 
of knot interactions in the high ejection rate cases (LJ0.5 and 
LJ2) and the absence in the low ejection rate case (LJ8) over the 
timescale considered. In the light of the above discussion, we 
suggest, therefore, that the ejection rate might be constrained if 
the frequency of mutual interactions among ejected blobs can 
be determined from observations. The constrain of the ejection 
rate could be also interesting in the context of the accretion sce- 
nario. In fact it is believed that the ejection/accretion phenom- 
ena are intimately related, but a full description of this relation 
is still lacking. 

5. Conclusions 

We investigated the physical mechanism determining the irreg- 
ular knotty structure of protostellar jets commonly observed in 
the optical band, and re cently also detected in X-rays (e.g. HH 
154. lFavata et al. 2006). To this end, we developed a numerical 
model describing a pulsed jet, consisting of a chain of blobs, 
each ejected along the jet axis with random initial velocity, and 
traveling in an initial homogeneous ambient medium. Our sim- 
ulations represent the first attempt to model a pulsed jet that 
simultaneously includes the thermal conduction and the radia- 
tive losses processes. In particular, we find that the thermal con- 
duction plays a crucial role in damping out the hydrodynamic 
instabilities that would develop within the cocoon. Since these 
instabilities contribute to break the jet, the thermal conduction 
makes the jet to propagate to larger distances. 

The jet is modeled both as an initial light jet (a scenario 
which seems likely to describe the o bservations of HH 154, 
both in X-rays and in the optical band: lBonito et al.|[2008|) and 
as an initial heavy jet. In all the cases, we find that the inter- 
actions among subsequent ejected blobs are a common feature 
to our simulations. As a result of these interactions, each new 
ejected blob travels through a strongly inhomogeneous medium 
characterized by density and thermal conditions largely varying 
along the jet axis; as a result, the density contrast between the 
blob and the surrounding medium varies during the jet/ambient 
evolution and the different pre- shock conditions leads to strong 
variations of the shock velocities and of the post-shock proper- 
ties. 

The mutual interactions among the blobs lead to a variety of 
density structures withi n the jet not descr ibed by current mod- 
els of jets (however see Yira k et al . 2009). These structures de- 
termine a rather complex morphology of the jet. The most rele- 
vant are: single high speed knots, showing a measurable proper 
motion in nice agreement with optical and X-rays observations; 
irregular chains of knots aligned along the jet axis, and possi- 
bly interacting with each other, consistent with observations of 
several protostellar jets; reverse shocks interacting with outgo- 
ing knots; oblique shock patterns produced by the reflection of 
shocks at the cocoon surrounding the jet. We also found that 



10 



R. Bonito et al.: Generation of radiative knots in a randomly pulsed protostellar jet 



the range of velocity of the knots formed from mutual interac- 
tions among bl obs is consistent with th e typical range of values 
observed (e.g. lReipurth & BallvlEoOll) . 

We explored the effects of the ejection rate of the blobs on 
the evolution of the density and temperature spatial distribu- 
tions of the jet/ambient interaction. The frequency of interac- 
tions between formed knots depends on the ejection rate. We 
suggest, therefore, that the comparison between our model re- 
sults and the observations may give some hints on the ejection 
rate of protostellar jets. Since the opposite and complementary 
processes of accretion of mass onto a star and ejection of mass 
from a star are commonly assumed to be related, from the con- 
straints on the ejection rate we may also infer important physi- 
cal information on the accretion rates. 

We expect some details of our model to depend on the as- 
sumption of the exponential distribution for the ejection veloc- 
ity of each blob (determined with a given random seed). From a 
quantitatively point of view, the single blob velocity can change 
with a different random seed and also the energetic and dy- 
namics of particular structures is expected to vary if different 
velocity distributions are chosen. Nevertheless, in any case we 
expect to obtain the main features discussed here, namely the 
formation of shock fronts, the collision among knots due to 
blobs ejected at different epochs with different speed, the inter- 
action among knots and reverse shocks. 

The model of a randomly ejected pulsed jet is, therefore, 
a promising scenario to explain the physical mechanism lead- 
ing to the irregular knotty structure observed within protostel- 
lar jets in a wide range of wavelengths, from optical to X-rays. 
Further development of this work will be reported in a forth- 
coming paper (Bonito et al., 2009, in preparation) where we 
will investigate the morphology and variability of X-ray emit- 
ting structures formed within the randomly ejected pulsed jet. 
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